'''
Created on Jul 19, 2013

@author: a.renduchintala
'''
from numpy import *
import math
# covariance matrix
sigma = matrix([[2.3, 0, 0, 0],
           [0, 1.5, 0, 0],
           [0, 0, 1.7, 0],
           [0, 0, 0, 2]
          ])

# mean vector
mu = array([2, 3, 8, 10])

# input
x = array([2, 3, 8, 10])

def norm_pdf_multivariate(x, mu, sigma):
  size = len(x)
  if size == len(mu) and (size, size) == sigma.shape:
    det = linalg.det(sigma)
    if det == 0:
        raise NameError("The covariance matrix can't be singular")

    norm_const = 1.0 / (math.pow((2 * pi), float(size) / 2) * math.pow(det, 1.0 / 2))
    x_mu = matrix(x - mu)
    inv = sigma.I        
    result = math.pow(math.e, -0.5 * (x_mu * inv * x_mu.T))
    return norm_const * result
  else:
    raise NameError("The dimensions of the input don't match")

print norm_pdf_multivariate(x, mu, sigma)


ss = matrix([[1, 0], [0, 1]])
mus = array([ 5.21943697 , 0.72891909])
xs = array([1, 2])

print norm_pdf_multivariate(xs, mus, ss)
'''
from scipy import stats
mynormnd = stats.norm(mu,sigma)
det = linalg.det(sigma)
size = len(x)
norm_const = 1.0/ ( math.pow((2*pi),float(size)/2) * math.pow(det,1.0/2) )
print prod(diag(mynormnd.pdf(x))) * norm_const
'''
